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Analytic formulas of Minkowski functionals in two-dimensional random fields are derived, including effects 
of second-order non-Gaussianity in the presence of both the bispectrum and trispectrum. The set of formulas 
provides a promising method to constrain the primordial non-Gaussianity of the universe by temperature fluctu- 
ations in the cosmic microwave background radiation. In a case of local-type non-Gaussianity, the Minkowski 
functionals are analytically given by powers of quadratic and cubic parameters, / NL and gr NL . Our formulas are 
not restricted to this particular model, and applicable to a wide class of non-Gaussian models. The analytic 
formulas are compared to numerical evaluations from non-Gaussian realizations of temperature maps, showing 
very good agreements. 
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I. INTRODUCTION 

Detecting the primordial non-Gaussianity will play a key 
role in discriminating the models of the early universe. All 
models of inflation predict the primordial non-Gaussianity to 
some extent in principle. While a simple inflationary model 
with a single slow-rolling scalar field produces too small non- 
Gaussianity to be observed 01-0], many other models such as 
the curvaton scenario 0-0], models with non-standard kinetic 
term ITolfTlll . cyclic or ekpyrotic universes without inflation 
lfl3tl . etc. predict large non-Gaussianity which could be de- 
tected in any future lfl4ll . 

The statistical properties of random Gaussian fields are 
completely characterized by the two-point correlation func- 
tion in configuration space, or the power spectrum in Fourier 
space. Thus the information of non-Gaussianity is con- 
tained in higher-order polyspectra, such as the bispectrum, 
trispectrum, and so forth. Recently, much work is devoted 
to analyzing the temperature fluctuations in the cosmic mi- 
crowave background (CMB) to constrain the primordial non- 
Gaussianity fl5l - l20ll . So far most of the constraints are placed 
to a phenomenological parameter /nl, which is largely re- 
sponsible for the primordial bispectrum. However, the bispec- 
trum is not sufficient to fully discriminate the wide variety of 
models of the early universe. Some models, such as the curva- 
ton scenario i2ll - l23ll . ekpyrotic universes ll24ll . etc. can gen- 
erate a large trispectrum with a relatively small bispectrum. 
The trispectrum is often characterized by phenomenological 
parameters #nl (and tnl in some models). When £/nl is large 
enough such that <?nl » f^ h , the bispectrum is not signifi- 
cant and the primordial non-Gaussianity can be detected only 
through the effect of trispectrum. 

Straightforward measurements of the polyspectra become 
progressively complicated for higher-order statistics beyond 
the bispectrum, since the higher-order polyspectra have 
many arguments with complicated dependence on shapes and 
scales. Even though direct calculations of the higher-order 
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polyspectra may be primary methods to constrain higher-order 
non-Gaussianity, it is desirable to have as many alternative 
methods as possible. Among others, the set of Minkowski 
functionals (MFs) I25l - f27ll has been proved to be a sim- 
ple and powerful tool in analyzing non-Gaussianity in CMB 

111 [MIS EMI. 

Since the hierarchy of polyspectra contains all the informa- 
tion of the statistical property of random field, the expectation 
values of MFs should be expressible by the polyspectra. Such 
relations were analytically derived for weakly non-Gaussian 
fields in the lowest-order approximation 13411 . It was shown 
that deviations of the non-Gaussian MFs from the Gaussian 
MFs are proportional to a linear combination of the bispec- 
trum. The analytic formulas are incorporated with the CMB 
bispectrum [35] and are applied to the analysis of CMB, suc- 
cessfully constraining non-Gaussianity parameter /nl without 
relying on massive numerical simulations to generate non- 
Gaussian CMB maps JMS. 

So far the analytic formulas of non-Gaussian MFs are 
known only in the lowest-order approximation, which we 
call first-order non-Gaussianity. The contributions from the 
trispectrum and other polyspectra to the geometrical de- 
scriptors, such as the MFs, appear in higher-order non- 
Gaussianities ll34l l39l |4QTI - To extract any information from 
MFs beyond the bispectrum, such as c/nl and tnl, it is par- 
ticularly useful to have analytic formulas with higher-order 
approximations. The purpose of this paper is to derive the 
analytic formulas with second-order approximation of non- 
Gaussianity, including the effects up to the trispectrum. These 
primary formulas are applicable to a wide class of non- 
Gaussian fields in 2-dimensional space. We also focus on 
how the formulas are applied to analyses of CMB temperature 
fluctuations. Assuming a local model of non-Gaussianity, we 
derive the explicit dependence of the MFs on parameters /nl 
and g NL (and r NL ). 

This paper is organized as follows: In §2, the analytic for- 
mulas of MFs with second-order non-Gaussianity are derived. 
Non-Gaussianity parameters to describe the analytic MFs are 
explicitly given by the bispectrum and trispectrum in §3. In 
§4, the polyspectra in the local model of non-Gaussianity are 
summarized for applications to the analytic MFs. Numerical 
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tests of the analytic formulas are given in §5, and our results 
are summarized in §6. 



II. ANALYTIC MINKOWSKI FUNCTIONALS WITH 
SECOND-ORDER NON-GAUSSIANITY 



The MFs are morphological descriptors to study statistical 
properties of random fields. For a smooth scalar field fix) 
with zero mean, (/) = 0, the excursion set over a given thresh- 
old v, Q — \x\f{x) > vcr} are defined, where cr = </ 2 }^ 2 is 
the rms of the field. The MFs for a set Q in a 2-dimensional 
space with smooth boundary dQ are given by B27I1 



Kdl, 



Vb(v) = f da, V,(v) = - f di, V 2 (v) = i f / 

Je 4 Jag 27r Jsg 

(1) 

where afa and d£ denote the surface element of Q and the line 
element along dQ, respectively, and k is the geodesic cur- 
vature on the boundary dQ. Apart from numerical factors, 
Minkowski functionals Vo, V\, V2 correspond to the area of 
Q, the length of dQ, and the curvature integral of dQ, respec- 
tively. According to the Gauss-Bonnet theorem, 2n times the 
Euler characteristic \ of Q is equal to the integral of geodesic 
curvature of dQ plus the integral of Gaussian curvature of Q. 
Therefore, when the excursion set is defined on the sphere of 
radius R, the Euler characteristic is given by a linear combina- 
tion of the Minkowski functionals via^ -V2 + Vo/2nR 2 B27I1 . 
In a flat space (R — > 00), % = y 2 . 

The expectation values of MFs in 2D are generically given 
by analytic forms 113411 



V k (v) = 



1 



0)2 



( 27r )(*+D/2 U2 _ kUk [^a- 



v k (v), (2) 



where cr\ = ((Vm) 2 ) 1 ^ 2 is the variance of the gradient field, 
and aj/t = n k ^ 2 /T(k/2 + 1) is the volume of the unit ball in k 
dimensions, i.e., cog = 1, a>i = 2, <±>2 = n. For a random Gaus- 
sian field, the MFs are analytically given by Tomita's formula 
[41], which corresponds to 



v k (v) = Hi-i(v), 



where 



H„(v) = e 



v 2 /2 



d 
"dv 



„-v 2 /2 



(3) 



(4) 



are the Hermite polynomials, and we define a function 

H_ x {y) ee e v2 ' 2 £ dve^' 2 = ^|erfc |-^| , (5) 

for k = in Eq. ©. 

Throughout this paper, we assume hierarchical orderings of 
the higher-order correlators, (f") c ~ 0(cr 2 "~ 2 ), where (• • • )c 
denotes the connected part, or the cumulant. In such a case, 
the reduced MFs Vk(v) are expanded by the variance cr as 



v k - v k +v k cr + v k cr + 



(6) 



where i»P = H^\(v) is the Gaussian contribution to the MFs. 
Primordial non-Gaussianities generated by most of the infla- 
tionary models satisfy the hierarchical orderings which we as- 
sume. In the previous work 13411 . the first-order corrections to 
the MFs in a non-Gaussian field are analytically derived: 



kSi 



k(k-l)Sn 



H k - 2 (v), (7) 



where 



S = 



(f)_ 

cr 4 ' 



Si = 



</ 2 v 2 /) 

cr 2 cr\ ' 



2<|V/| 2 V 2 /> 



5 n = 



(8) 



are the skewness parameter and its derivatives. 

Fundamental techniques to evaluate the next-order term, 
vf\ are mostly found in j34ll . We outline below the higher- 
order extensions of the calculation in 03411 . First we define the 
normalized cumulants of the field derivatives, 



M (n) 



2-2« 



(9) 



>(«) 



(which corresponds to the quantity denoted by M^...^ „, 
Il34l0 . where denotes all the field derivatives, (f^) - 
(f,f;uf;2,.f;U,f;22,f;i2,--\ and ";" indicates the covariant 
derivative. In the following we adopt a notation, /o = /, 
fi - /;i> fz = f;2> f\\ = /;ii> an d the corresponding index 
runs over fi = 0, 1, 2, 1 1. The Eq. © is of zero-th order in cr 
for hierarchical orderings. For any statistic which is locally 
defined by some function F of field derivatives f^, the expec- 
tation value of the statistic is expanded as (Eq. (22) of [34]) 



(F) = F G + - V M (3) F G cr 



- y m (4) f g 
+— y M^i „ m, ( , 3 ;, „ f, 
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(10) 



d, = 



d/df„ 



where F G = (F) c , F® m ... = (d^ ■ ■ ■ F) c 
and (• • • )g denotes the expectation value for the Gaussian field 
which has the same two-point functions as the non-Gaussian 
field we consider (see 13411 for detail). The first two terms on 
RHS of Eq. ( fTOb correspond to the Gaussian contribution and 
first-order non-Gaussian corrections, respectively. 
For the MFs, the function F is given by 1 3 Z 



F = 



'©(« - v), (k = 0), 

f5(n-v)|ui|, (k=l), 
— \5(u - v)6(u\)\u2\u\\, (k = 2), 



(11) 



where 0, 6 are the step function and the delta function, re- 
spectively, and = f pier. The quantities F G are explicitly 
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calculated to give 1 3 



Eq. ( fTUI ). we obtain the second-order corrections to the MFs: 



k—li—lz—lm 



Hn-i(y)5 w 8i 2 o8 m o, (k = 0), 

x ■ ^H„(v)hi^ 2 5t 2 oS m o, (k = 1), 

hh h h-2 [H„+\(v)5 m0 - H n (v)6 m i] (k = 2), 

(12) 

where q = <T\I y/2cr, and n,l\,l 2 ,m are numbers of 0, 1, 2, 

11, respectively, in a set of indices fii,fi2, The factor hi is 

defined by 



^ 2) (v) 



vf\v) 



4 2) (v) 



^(v) + |h 3 (v), 
1 



m 



(17) 



(18) 



J' Hl(v) + Ez^li H5(v) 

72 24 



(19) 



0, 



(/ : odd), 



I (-2) //2 7r- 1/2 r[(Z + l)/2], (/ : even). 



(13) 



For example, h-2 - 1, ft-i = 0, fto = 1, /?i = 0, h% = -1, 
/13 = 0, /14 = 3, etc. 

The factor M^^ 3 is expressible by the skewness and its 
derivatives of Eq. ® from rotational symmetry lf34tl . The non- 
zero components are 
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000 



M 
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and their permutations. For other combinations of the indices, 
^K"iiu2j"3 = 0. Substituting Eq. ( fT2l i and Eq. (IT4-b into the second 
term in the RHS of Eq. ([Tol l, the first-order corrections to the 
MFs of Eq. ©follow. 

It is straightforward to extend the above calculation to de- 
rive second-order corrections. We only need to calculate 4- 
point cumulants M^imu^u, ■ With similar considerations of ro- 
tational symmetry in the case of 3-point cumulants, we obtain 



M (4) = K M (4) 

JKI 0000 m 000(11) 
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M 022(ii) -1 (Kn-Km), 
< = 9 4 K m . 
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(ID - -9 K m, 
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M, 
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2222 



3q 4 K m , 



(15) 



For other combinations of the indices which cannot be ob- 



(4) 



tained by permutations in the above equations, 
Kurtosis parameter and its derivatives are defined by 



= 0. 



K = 



</ 4 >c 



K, = 



(fv7>c 

crV? 



2</|V/| 2 VY)c + <|V/I 4 )c 
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(T L (J X 



in 



<|V/I 4 )c 

2a- 2 o-\ 



Substituting Eq. (O-CB]) into the third term in the RHS of 



These equations are our primary results of this paper. In 
Eqs. ©, ©, 0, <fT7V([r9l. the MFs of non-Gaussian fields 
are analytically given by the skewness, kurtosis, and their 
derivatives in second-order approximations. The above for- 
mulas are applicable to any 2D random field which has hier- 
archical orderings of the higher-order cumulants. 



III. RELATION TO THE BISPECTRUM AND 
TRISPECTRUM 

Once the skewness, kurtosis and their derivatives are cal- 
culated from a cosmological model, the prediction of MFs 
follows from our formulas. The skewness and its derivatives 
are directly given by the bispectrum, and the kurtosis and its 
derivatives are directly given by the trispectrum. 

For a field on a 2D flat space, such as the flat-sky approx- 
imation of the CMB temperature fluctuations AT(9)/T, one 
can adopt the Fourier transform: 



T V ' J (In) 1 V ' 



(20) 



The Fourier coefficients of the smoothed field / is given by 
/(Z) = a(l)W(l), where W(l) is the window function of the 
smoothing kernel, where we assume spherically symmetric 
smoothing window. For a Gaussian window with smoothing 
angle # s , W(l) = e s ^ 2 . The power spectrum C, the bis- 
pectrum B, and the trispectrum T of the 2D field AT/T are 
defined by 

{a{h)a{l 2 )\ = (2n) 2 6 2 (h + l 2 )C(h), (21) 
(a(h)a(l2)a(l 3 )) c = (2tt) 2 < 5 2 (Z 1 + l 2 + h)B(h, l 2 , Z 3 ), (22) 
(a(li)a(l2)a(l 3 )a(l 4 )) c = {2n) 2 6 2 (l x + Z 2 + Z 3 + h) 

XT(h,l 2 ,h,l 4 ;l 12 ,l23), (23) 



(16) where 1/ = |Z,-|, Zy = |Z; + lj\. The translational symmetry 
is guaranteed by delta functions, and the rotational symmetry 
is taken into account in above definitions of arguments. By 
straightforward calculations using symmetric permutations in 
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integrals, one can show that 



K A = 



= [LA 
^1" 



1 A/ 



Z] rfZj Z 2 d/ 2 <#?12 p , , 
"Tj = &A(!ul2) 

In In In 

xB(h,l2,h2)W(h)W(h)W(h 2 ), 

h dh l 2 dl 2 h dh dO n dQn ~ 

KA{l\,h, hi) 



2n 2n 2n 2n 2n 
x T(l u h, h, h\ hi, li3)W(h)WQi)W(h)W(l4) 



where 



hi = {l\ + + 2hh cos On) 1 ' 2 , 
l 23 = (l\ + + 2hh cos fe) 1/2 . 
U = [l 2 12 + l W + 2l ih cos(6» 12 + 23 )] 1/2 , 



(24) 

(25) 
(26) 

(27) 
(28) 
(29) 



and tr = <r, (S A ) = (S,Sj,S n ), (S A ) = (§,§i,Su), (K A ) = 
(K, K\, Kn, Km), (Ka) = (K, K\, K\\, Km), 



K=l, Ki=- 



ft 
h 

'2q 2 
f 

__L 

2q 2 



Su ■ 
K a 



l\ (lj-2ll) 



Aq 4 



I 4 -4l 2 l 2 
hi ^'i'3 



16q 4 



Km — 



l\ 1+ M\(l\-l 2 n ) 
32q 4 



(30) 
(31) 
(32) 



The integration ranges are 2 < I, < oo and < 6jj < 2n when 
the monopole and dipole components are subtracted from the 
map [see Eqs. (l53T>-d55l> and explanations below]. 

For a field on a 2D sphere, such as the all-sky CMB map 
AT(6,(f>)/T, one can adopt the expansion by spherical har- 
monics: 



AT 



(9,$) = Yj a '™ Y T( e > 



(33) 



The harmonic coefficients of the smoothed field / is given by 
fim = cii, n Wi, where Wi is the window function. For a Gaus- 
sian window with smoothing angle 9 S , Wi = e~ /( ' +1)fl s/ 2 . When 
monopole and dipole components are subtracted as usually in 
the case of CMB map, Wo = W\ = 0. The power spectrum, 
bispectrum and trispectrum are defined by J4J-SJ1 



(a/ imi a/ 2 m 2 ) c - (-ly'fy/jdmj ,-m 2 Qi 
h h h 



(a hmi a hm a hm ) c = I ^ ^ £\B hhk 



(34) 
(35) 



= 2>iW 

1 M \ 



(h h 


M 




U 




\m\ m 2 


-M J 




n?4 


MJ 



where ui,JjJj 3 J is the Wigner 3y'-symbol, and the rotational 
symmetry is taken into account in the above definitions. The 



symmetries of the trispectrum T^J^L) are enforced by the re- 
duced trispectrum T}b(L) which is an arbitrary function of its 
arguments except that it must be symmetric against exchange 
of its upper and lower indices T^(L) = T^(L) Q. The 
construction is 



, (37) 



where {^j^} is the Wigner 6 /'-symbol, and 



+ (-l) h+u+L T[£(L) + (-i)'i + '=+'3+'47-^i (L). (38) 



h+k+h+kn-hh 



The following properties of spherical harmonics B45H are 
useful for our purpose: 



v 2 y ; m (Q) = + i)if(n), 

J" dDY^ (Ci)Y^(0) = (-\) mi S hk 8 mu . m2 , 

y™'(Q)y;: 2 (Q) = Yj(-iri hkh 



(39) 
(40) 



I3 ,lfli 



h h h |^3 (fi) (41) 

mi nil ni'i 



where 



(2/j + 1X2/2 + 1X2/3 + 1) //! l 2 M 

hhh = \l ^ ooo- (42) 



This quantity 7/,;,/ 3 is completely symmetric with respect to 
permutations of its arguments, non-zero only when h,l 2 ,h 
satisfy the triangular inequality and h + h + ^3 is an even 
number B45I1 . With repeated use of Eqs. <f39])— (f4Tb. one can 
calculate the required parameters. Introducing a notation, 



{/} = Kl + 1), 



(43) 



the results are 



(44) 
(45) 



S A = - — 7 V I1M3S AhhkBhhkWhWizWis, 

(46) 



h,h,hJ4,L 
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where 



{/l} + {/ 2 } + {/3 
6q 2 



(47) 



\h) 2 + \h) 2 + {hf ~ 2{h){l 2 ) - 2{h)\h) - 2\h)\h 



I2q 4 

{/l} + {/ 2 } + {/ 3 } + {/4} 



Phhm W 2 -W l} + fe})({/3} + fa}) 



Wq 4 

({h} + {l 2 } - {L}) ({h) + {l 4 } - {L}) 
32q 4 



(48) 
(49) 

(50) 

(51) 



The above forms of skewness parameters [Eq. ( l47l i. (l48t l were 
already appeared in 13511 . 

Resemblances of the above results to those of the flat space 
are obvious if the integrands of Eqs. d3Qb — (f32T> are sym- 
metrized with respect to l\,h, h, and U [conversely, one can 
desymmetrize the Eqs. d47ll— (fSTb to have the similar form with 
Eqs. d30b-d32bl. Noting the all-sky and flat-sky correspon- 
dence 146, H^, it is a straightforward exercise to show 
that the above all-sky equations reduce to those of flat-sky in 
the large-/ limit. Following ll44ll46ll47ll . but applying an im- 
proved approximation 



Kf(0,0) *(-!)" 



21+ 1 
4n 



Jimp 



(52) 



for 8 <k 1, I » 1, the correspondences between all-sky and 
flat-sky spectra are derived as 



B hhh ~ hhhB^l\ + -,l 2 + -,h + - 



(53) 
(54) 
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xP\h + 2^2+ 2 ;/3 + 2' Z4+ 2' L+ 2 



and 



T(lu h, h, k; In, I73) = P(h , h; h, U\ hi) 

+ P(h,h; h, U; in) + P(h,k; h, h; ha), (56) 

where l n = [I 2 + l\ + l\ + l\- l 2 2 - l\^ 12 . 

Since the all-sky multipole I and the flat-sky wavenumber 
|Z| are related by |Z| — £ + 1 /2, contributions of the multipole I 
to the all-sky summation is approximately represented by the 
flat-sky integration over the range i- 1/2 < |/| — 1/2 < (+1/2, 
i.e., € < \l\ < £+ 1. Thus, all-sky summations over^ = 2,3, . . . 
correspond to the flat-sky integrations with the limit \l\ > 2, 
as noted above. We confirm that the flat-sky approximations 
of Eqs. (l24b-(|26T) with the above correspondences numeri- 
cally reproduce the values calculated from all-sky formula of 
Eqs. d44li-(l4"6jl within several percent for 9 S < 100'. 



For numerical evaluations of the kurtosis and its derivatives 
by the summation of Eq. d46l ), the number of terms to add is of 
order 0(l^ ax ), where l max is the maximum multipole required 
for a given smoothing scale, e.g., l m . dx ~ several X The 
computational cost becomes progressively high for large l mdx , 
if the summation is naively performed. Efficient evaluations 
are necessary when the smoothing angle 9 S is small. In Ap- 
pendix|Aj numerical schemes for the efficient evaluations are 
summarized. 



IV. A SIMPLE EXAMPLE: THE LOCAL MODEL OF 
NON-GAUSSIANITY 

The analytic MFs are evaluated once the power spectrum, 
bispectrum and trispectrum are given. These spectra depend 
on models of primordial density fluctuations. As a simple ex- 
ample, we consider below the local model of non-Gaussianity, 
although our formulas are not restricted to this particular 
model. 

In the local model, the primordial curvature perturbations 
during the matter era is assumed to take the form rf42ll43ll48l — 

m 

®(sc) = 4>{x) + / NL (cf> 2 (x) - <0 2 >) + g Nh cf> 3 (x), (57) 

in configuration space, where <p is an auxiliary random Gaus- 
sian field. The comoving curvature perturbation £ is given by 
f = 3<D>/5. The CMB fluctuations generated by the curvature 
perturbations have the harmonic coefficients 



aim = 4n(-i) 



-(TV 



(58) 



where <J)(fc) is the Fourier transform of the primordial cur- 
vature perturbation <t>(x), and gP(k) is the radiation transfer 
function. 

The bispectrum and trispectrum of CMB in the local model 
of Eq. (|57]) are derived in literature [SHH]: 



(55) Bi t i 2 i 3 - 2/ NL /;, 



s 



r dra h (r)fi h (r)P h (r) + eye. 



(59) 



rfdn r\dr 2 F L {ri , r 2 )a h (r\ )/3 h (n )a h (r 2 )p u (r 2 ) 



+ g NL j r 2 dr/3, 2 {r)j3 u (r) [a h (r)0 h (r) + /3 h (r)a h (r)] 



where 



F L (n,ri) = 4n 
a,(r) 



fair) = 4n 



I 

I 



k 2 dk 
2^" 

k 2 dk 



P0(k)j L (k ri )j L (kr 2 ), 



2n 

k 2 dk_ 
2^2 



2 g™(k)M*r), 



PWgP{k)Mkr), 



(60) 

(61) 
(62) 
(63) 
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and P$(k) oc k"< 4 is the primordial power spectrum of 0. In calculation: 
some extended models, the factor 4/^ L in the trispectrum is 



replaced by 4/|Jj L — > 25tnl/9, and tnl is considered as an 
independent parameter 115111 . 

The variance parameters cr and en are also affected by non- 
Gaussianity. Therefore we need to evaluate the corrections 
to the power spectrum. The non-Gaussian corrections to the 
parameters cr, cr\ mostly affect the normalization of MFs. In 
practice, the normalization suffers contamination from obser- 
vational effects such as pixelization and/or boundary effects, 
and therefore is not usually used for extracting cosmological 
information. Nevertheless, we include the effect for theoret- 
ical consistency here. By substituting Eq. d58"b into Eq. d34l i 
for a local model of Eq. ( fSTt , performing angular integrations, 
we obtain 



Ci 



where 



1 + 6g N L r dr^(r)\Q 



2/nl j ^ dry if) [^(r)f , (64) 

(66) 



' k 2 dk . 



P^k)j (kr), 



7i(r) = (An) 1 J ^ [g ( P(k)f j (kr). (67) 



The quantities C;, ^(r) are the angular power spectrum from 
the Gaussian component, and the spatial correlation function 
of 0, respectively. Substituting Eqs. d59l l, ( |60i >, ( |64l into 
Eqs. d44l>-(l4"6i> with Eqs. (1371 ). (l38l l. all the parameters nec- 
essary to calculate the analytic MFs are evaluated for a model 
of local non-Gaussianity with a full transfer function g^\k). 

For the Sachs-Wolfe effect HH, which is valid only for 
small multipole moments (I <c 100), the radiation transfer 
function is simply given by g^(k) — ji(kr x )/3, where r, is 
the comoving distance to the last scattering surface. In this 
limit, the reduced bispectrum and trispectrum have the forms 

Emm 



(68) 



a hhh = -6/NL [C h t h +L, 2 L h +L h L h )h\hh. 



7^ (L) = 9C sw c sw 



where 



x f^ L Cf w + g Nh (C™ + Cf 3 w )] l hhL l hUL , (69) 



Cf w = ^ J k 2 dkP^k) [jtikr*)] 2 , (70) 

is the power spectrum of the Sachs-Wolfe effect from the 
Gaussian component. In the same limit, non-Gaussian cor- 
rections to the power spectrum are also obtained after some 



1 + 270NLy (2L+1)c S W 



In 



C 



svv 



21 



hhl 1 



r SW r SW 



(71) 



V. COMPARISONS WITH NUMERICAL MINKOWSKI 
FUNCTIONALS 

To illustrate how well our formulas work, we compare the 
analytic MFs with numerical MFs from realizations of the 
non-Gaussian CMB map. Since the purpose of this com- 
parison is to show the performance of analytic formulas, we 
simply consider the non-Gaussian map in the Sachs-Wolfe 
limit. Implementations with a full radiation transfer func- 
tion will be used in future applications to real data. In the 
Sachs -Wolfe limit, the temperature fluctuations are given by 
AT sw (0,(p)/T = -<$(r t ,0,(p)/3, and thus generated by the 
nonlinear mapping of Gaussian fluctuations ATc/T at each 
position in the sky: 



AT; 



sw 



AT G (AT G 
— -Vnl 



AT C 



- (monopole + dipole), (72) 

where monopole and dipole components are subtracted from 
the resulting map. We use the HEALPix package 15 311 to 
generate the non-Gaussian map, apply the smoothing win- 
dow function, and calculate the first- and second-derivatives 
of the temperature field. The MFs are evaluated from the field 
derivatives using the numerical algorithm given in l27l l35ll . 
We assume a simple Sachs-Wolfe power spectrum of 1(1 + 
l)Cf w /27r = 10~ 10 for I < 128 and Cf w = 0, otherwise. 
Without this cut-off, the resulting power spectrum would log- 
arithmically diverge in Eq. dTTl i (the actual power spectrum 
has natural damping in high-/ regime). We adopt the Gaussian 
filter W/ = e -'( /+1 ^/ 2 with the smoothing radius of 6 S = 100', 
and non-Gaussianity parameters /nl = 10 2 and <7nl = 10 6 . 
We generate 100,000 realizations of the non-Gaussian map, 
numerically calculate MFs in each realization, and finally av- 
erage over the realizations. 

In Fig. [1] the analytic and numerical MFs are compared. 
MFs as functions of the threshold v are plotted in upper pan- 
els. As noted in [35], precise values of overall normalizations 
are affected by numerical artifacts such as pixelization, bound- 
ary effects etc. even in the Gaussian random field at the sub- 
percent level. However, deviations from the Gaussian shape 
of MFs as functions of the threshold v are in very good agree- 
ment between numerical and analytic MFs. In lower panels, 
differences between non-Gaussian MFs and Gaussian MFs, 
divided by the maximum amplitudes, are compared. For nu- 
merical data, MFs from Gaussian realizations are calculated 
and subtracted from those of non-Gaussian realizations. We 
find the agreements are extremely well. Note that there is no 
fitting parameter at all in the comparisons. 
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FIG. 1: Comparisons between analytic MFs (solid lines) and numer- 
ical MFs (symbols) in the Sachs-Wolfe limit. We use 9 S = 100' and 
/nl = 10 2 , gNL = 10 6 . Contributions from /nl (dotted lines) and g^L 
(dot-dashed lines) are also shown. Total MFs are given by sums of 
the two contributions in such a weak non-Gaussianity regime. 



Contributions from the parameters /nl and g^L are shown 
in thin dotted lines and dot-dashed lines, respectively. Phases 
of oscillating patterns are different from one another. This 
is because the higher-order non-Gaussianity involves higher- 
order Hermite polynomials as given in Eqs. (0 and (TT7T) — <fT9b- 
Since the Hermite polynomials are orthogonal functions, the 
effect of each parameter can be distinguished with such obser- 
vations. 



VI. SUMMARY 

In this paper, the analytic MFs with second-order non- 
Gaussianity, including effects from the bispectrum and 
trispectrum, have been derived. As long as the higher-order 
cumulants obey hierarchical structure, (f) c ~ 0(<t 2 "~ 2 ), de- 
viations from the Gaussian predictions of MFs are expanded 
by cr. First several cumulants are relevant to the non-Gaussian 
corrections. Up to the second-order in cr, the non-Gaussian 
corrections are expressed by the skewness, kurtosis, and their 
derivatives, and thus by the bispectrum and trispectrum. These 
fundamental results are general and applicable to any 2D ran- 
dom field which has hierarchical orderings of higher-order cu- 
mulants. In cosmology, applications to the CMB temperature 
fluctuations are of great importance. Once the form of bispec- 
trum and trispectrum in CMB are given by a model of the early 
universe, the analytic MFs are evaluated by our formulas. We 
have illustrated how the analytic MFs of CMB are calculated, 
employing the local model of non-Gaussianity, as an example. 
We have compared the analytical MFs with numerical calcula- 
tions of non-Gaussian MFs, and found very good agreements. 
Applications to real data of CMB map are now in progress and 
will be reported in future work. 
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Appendix A: Evaluations of Skewness, Kurtosis and Their 
Derivatives 

When the Eq. d46*l i is straightforwardly evaluated, the num- 
ber of terms to add is of order <9(/ max ), where Z max is the maxi- 
mum multipole moment required for a given smoothing scale, 
e-g-, 'max ~ several x 6~ l . For parameters K and K\, summa- 
tions over L for terms with 6j-symbols can be analytically 
performed, but we still need numerical summations of 6j- 
symbols for parameters Ka and Km- The computational cost 
is progressively high for large Z max , and therefore efficient cal- 
culations are necessary when the smoothing angle is small. 

First of all, one should take advantages of symmetries and 
constraints to save the number of additions. For the skewness 
parameters, the number of addition is reduced by a factor of 
about six, using the substitution 



=h h 



h,h,h 



=l 2 <h 



• 3 Z +6 L ' (A1) 

h<h=h h<h<h 
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in Eq. d45l l. The summation is taken only when l\ + I2 + h is 
an even number and I1J2J3 satisfy the triangular inequality, 
h—l\<h< l\ +h because of the factor //,/,/,. For the kurtosis 
parameters, the number of addition is reduced by a factor of 
about eight, using the substitution 

Z = Z +2 Z +2 X +2 z 

h,h,h,U h=i 2 ;,=;, ;,=; 2 /,</, 
/ 3 = {4 '3 = '4 '3 < U h = U 
h =h h< h h =h h = «3 

+ 4 Z +4 Z +4 Z +8 Z ' < A2 > 

l\ = l 2 l\ < h h < h h < h 

/ 3 < l 4 l 3 = I 4 Z 3 < / 4 / 3 < [4 

l\ < 1$ l\ < l 3 l\ = / 3 l\ < / 3 

in Eq. (l46l l. The summation is taken only when l\ +h + h + U 
is an even number and I4 < /1+Z2+/3 is satisfied because of the 
factor h^hh-iinh- The index L in Eq. d46b runs over every other 
integers in the range [max(|/i — fcl, I/3 — Z4I), min(Zi +h,h + W] 
with the same parity as that of l\ + 1%. 

Secondly, since the trispectrum is not a strongly oscillating 
function in usual cases, one can sparsely sample the multi- 
poles in the sum for large values of l\, . . . ,U. Acoustic oscil- 
lations in the CMB are sufficiently mild in this respect. There 
is a caveat that the summation over L' in Eq. d37l > should not be 
sparsely sampled, since the 6 j-symbol is a strongly oscillating 
function with respect to its arguments. Sparse samplings for 



large Z,'s with appropriate weights save enormous amounts of 
time for the calculation. 

Instead of calculating the exact summation, one can use 
the flat-sky approximation of Eq. (1261 . and evaluate the five- 
dimensional integral by, e.g., the Monte-Carlo algorithm. 
Once the reduced trispectrum T'^iL) is given, Eqs. d38l , d551 l. 
(l56i > give the trispectrum T(li,l2, h, U\ In, hi) in the flat-sky 
approximation. There is no need for summing 6j'-symbols 
in this procedure. The integration ranges of Eq. d26l i are 
2 < h,h,h,h < 00 , < #12, #23 - 2n when the monopole 
and dipole components are subtracted from the map. This ap- 
proximation is valid when most contributions to the kurtosis 
parameters come from the flat-sky regime (I » 1). In realis- 
tic spectra, this condition is satisfied because multipoles near 
the acoustic peak at I ~ 200 dominantly contribute. Strictly 
speaking, the flat-sky integration is valid when all fs are large. 
Therefore, for a more precise approximation, one can use the 
all-sky summation of Eq. (l46*T l only when, e.g., l\ < 20 [note 
that l\ < I2, h,U m Eq. (1A2H . and otherwise adopt the flat-sky 
integration ofEq. d26ll over the range l\, h, h, U > 21. 

Unless sub-percent level accuracies are required, experi- 
ences show that the simplest flat-sky integrals over all ranges 
of multipole with appropriate choice of integration limits suf- 
fice for evaluations of skewness and kurtosis parameters in 
realistic spectra. 
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